International Journal of Electrical and Computer Engineering (IJECE) 
Vol. 12, No. 3, June 2022, pp. 2553~2563 
ISSN: 2088-8708, DOI: 10.1159 I/jece.v12i3.pp2553-2563 O 2553 


Brain tumor visualization for magnetic resonance images using 
modified shape-based interpolation method 


Dina Mohammed Sherif El-Torky!, Mohamed Ismail Roushdy?, Maryam Nabil Al-Berry', 
Mohammed Abd El-Megeed Salem? 


'Faculty of Computer and Information Sciences, Ain Shams University, Cairo, Egypt 
Faculty of Computers and Information Technology, Future University in Egypt, Cairo, Egypt 
*Faculty of Media Engineering and Technology, German University in Cairo, Cairo, Egypt 


Article Info ABSTRACT 

Article history: 3D visualization plays an essential role in medical diagnosis and setting 
. treatment plans especially for brain cancer. There have been many attempts 

Received Dec 18, 2020 for brain tumor reconstruction and visualization using various techniques. 

Revised Dec 28, 2021 However, this problem is still considered unsolved as more accurate results 

Accepted Jan 10, 2022 are needed in this critical field. In this paper, a sequence of 2D slices of 


brain magnetic resonance Images was used to reconstruct a 3D model for the 
brain tumor. The images were automatically segmented using wavelet 
Keywords: multi-resolution expectation maximization algorithm. Then, the inter-slice 
gaps were interpolated using the proposed modified shape-based 
interpolation method. The method involves three main steps; transferring the 
: : binary tumor images to distance images using a suitable distance function, 
Shape-based interpolation interpolating the distance images using cubic spline interpolation and 
Visualization thresholding the interpolated values to get the reconstructed slices. The final 
tumor is then visualized as a 3D isosurface. We evaluated the proposed 
method by removing an original slice from the input images and 
interpolating it, the results outperform the original shape-based interpolation 
method by an average of 3% reaching 99% of accuracy for some slice 
images. 
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1. INTRODUCTION 

Brain cancer is one of the most life-threatening diseases to human beings. More than 22,000 people 
in the United States of America are diagnosed with brain tumors every year [1]. The treatment process 
involves one or more of surgical resection, radiotherapy, and chemotherapy. For the oncologist to decide a 
suitable treatment plan, there must be an accurate and reliable diagnosis of the pathology, location, shape, 
and size of the tumor. The diagnostic techniques are mainly divided into invasive techniques, e.g., biopsy, 
and non-invasive techniques such as medical imaging tools. 

Magnetic resonance imaging (MRI) is an important imaging tool for medical diagnosis especially 
for lesions. The images can show the lesion properties when injecting the patient with Gadolinium-based 
contrast agent [2]. However, the MRI yields a set of cross-sectional 2D images that need further 
interpretation from an expert radiologist to construct a 3D model of the tumor from these slices. This can 
result in an inaccurate diagnosis because the analysis depends on the experience of the radiologist [3]-[5]. 
Moreover, to avoid the slices from interference with each other (crosstalk artifacts), there must be an inter- 
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slice spacing during the image acquisition [6] which causes loss of details in the gap spaces and extra effort 
from the radiologist to analyse the images. 

As a result, many attempts are being made to provide an automated visualization for brain tumors. 
We propose a modified shape-based interpolation technique to estimate the missing information in the gap 
slices then reconstruct a 3D model of the tumor that can be easily viewed and interpreted. First, the tumor is 
segmented from each MRI slice. Then the tumor images are converted to distance images which are used to 
interpolate the in-between missing slices. Finally, the tumor volume is rendered and displayed in a 3-D view. 

The rest of the paper is organized; section 2 reviews the most popular methods that are currently 
used for tumor 3D reconstruction. Section 3 explains the methods used for tumor segmentation and 
extraction. The proposed method is presented in section 4. Results are shown in section 5 and the conclusions 
are demonstrated in section 6. 


2. RELATED WORK 
3D reconstruction is implemented using various methods that are based on diverse concepts. In this 
section, we will discuss the most widely-used methods [3]. 


2.1. Delaunay and alpha-shapes 

Delaunay methods are mainly based on extracting tetrahedron surfaces from an initial point cloud. 
The concept of 3D reconstruction using alpha-shapes was first introduced by Edelsbrunner and Mucke [7]. 
Using a finite set of points S and the parameter alpha, the resulting alpha-shape of S is a polytope, which is a 
generalization in any number of dimensions, of the three-dimensional polyhedron that is not necessarily 
convex or connected. If the alpha value is large, the alpha shape of S is similar to its convex hull. By 
decreasing the value of alpha, the non-convex details of the shape are formed. Edelsbrunner and Mucke [7] 
modified the technique by forming a surrounding sphere smaller than alpha then removing all tetrahedrons 
that lie inside it. External triangles of the remaining tetrahedron form the final surface. Al-Tamimi et al. [8] 
used alpha-shape theory for 3D reconstruction of brain tumors using slices of MR images. Most 
delaunay-based methods have the advantage of accurately fitting the surface defined by the original point 
cloud. On the other hand, their performance is strongly affected with noise. 


2.2. Zero-set methods 

Zero-set methods are also known as implicit reconstruction methods. They use a distance function 
for surface reconstruction. A polygonal representation of the object is formed by extracting a zero-set using a 
contour algorithm. Therefore, surface reconstruction from a point cloud is simplified to computing the 
suitable function f which equals zero for the sampled points and does not equal zero for the rest of the 
points [3]. The most popular zero-set method is the marching cubes algorithm which was proposed by 
Lorensen and Cline [9]. Guo et al. [10] introduced an improved marching cube algorithm by combining the 
seeded region growing and the standard Monte Carlo (MC) algorithm. Many other methods were proposed 
based on Zero-set methods such as Kazhdan and Hoppe [11] and Walder et al. [12] The drawback of these 
methods is the loss of geometrical accuracy in high-curvature areas such as corners and edges of the re- 
constructed shape. 


2.3. Point-based methods 

From the popular point-based methods is the voxel-grid filtering. It is based on the idea of reducing 
the number of points by sampling the input space using a grid of 3D voxels then a centroid is selected for 
each voxel to represent all the points [13], [14]. Angelopoulou et al. [15] used growing neural models to 
automatically landmark the target volume sections and construct a 3D model. The shortcoming of these 
methods is that it is impossible to determine the final number of points representing the surface. Methods that 
depend on deep learning and neural networks as in [16] and [15] were also proposed. Although that these 
methods produce highly efficient output, they have a drawback of finding the suitable training data especially 
with brain tumors that have various shapes, sizes, and consistencies. 


2.4. Interpolation-based methods 

To avoid interfering of slices during MRI acquisition (crosstalk artifact), there must be a gap 
between the images. The interpolation-based methods idea is based on calculating the image values in the 
gap area to reconstruct the 3D model accurately. These methods can be divided into scene-based and 
object-based methods. The simplest scene-based method is the linear interpolation. More methods were then 
proposed such as cubic spline and other polynomials in interpolation of medical images [17]. The kriging 
method was also used to interpolate grey values of medical images [18]. The scene-based methods apply 
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intensity averaging on the neighboring slices without considering the shape feature distortion. This causes 
blurriness to the object boundaries of the resulting interpolated slices. 

Several object-based interpolation methods were proposed to avoid the shortcomings of the 
scene-based methods. They considered object information in the given input slices to guide the interpolation 
procedure. An efficient type of object-based methods is the shape-based interpolation [19], this method uses 
distance transform functions to build distance images from the given slices. Then, instead of interpolating the 
intensity values, it interpolates the displacement preserving the geometric changes of the objects more 
accurately. The shape-based methods have become widely used due to their efficiency. Yet, they have the 
problem of disability to effectively deal with object having holes, large offsets, or heavy invagination. 

In [20], multi-resolution registration was used for image interpolation. Farias et al. [21] used 
interpolation techniques to generate intermediate slices of medical images. Ting ef al. [2] applied surface 
reconstruction using Hermite surface interpolation for breast cancer. Our proposed method uses shape-based 
interpolation to estimate the intermediate gap slices of the brain tumor volume. 


3. PROPOSED METHOD 

Our system is based on two main processes: tumor segmentation of brain MR images and 
intermediate slices interpolation of binary tumor images. Various segmentation methods have been proposed 
throughout the past years. These methods could be very simple such as Thresholding, or region-based as 
region growing and watershed, or supervised classification as in artificial neural networks, or statistical-based 
such as expectation maximization and its modifications [22]-[26]. We used wavelet multi-resolution 
expectation maximization algorithm (WMEM) [27] proposed by Salem [28], [29] for tumor segmentation. 
Besides being an unsupervised method, WMEM showed very accurate segmentation results for the brain MR 
images. Figure | shows the process flow diagram of the proposed method. 


Binary Tumor 
images 


Tumor Distance 


Segmentation 


Extraction Transform 


vryvyyy 
sh eh A Be 


Brain MRI series | Distance images 
AAAAS 


Gap Slices 
Interpolation 


Interpolated 
vrrvyy Interpolated Distance images 
Binary images AAAAAAAAAAA 


7 


Interpolated 


Slices 
Binarization 


Visualization 


fwevevere 


3-D Model 3-D Reconstruction 


Figure 1. Process flow diagram of the proposed method 


3.1. Multiresolution-based segmentation 

The WMEM is basically a modification of the classic expectation maximization (EM) algorithm. 
The EM consists of two major steps: expectation step (E-step) and maximization step (M-step). These two 
steps are repeated until convergence, i.e., the difference between the parameters in the previous and the 
current iterations decreases until it reaches a given threshold. The EM then uses a classifier that assigns a 
class membership to a pixel i depending on its intensity x;. The class to which the pixel is assigned is the one 
having a parameter vector that maximizes the Gaussian density function. The WMEM uses a Haar wavelet 
transform function that produces four outputs: image approximation, and vertical, horizontal and diagonal 
details of the image [30]-[32]. Two levels of approximations were produced to form a parent and a 
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grandparent of the original image. Since the edge details are lost during formation of the parent and 
grandparent images, these pixels are usually misclassified. 

To overcome this, a mask containing all the edges of the image is formed by adding the vertical, 
horizontal and diagonal details images resulted from the wavelet analysis. This mask is used to exclude the 
pixels laying on the edges of the image from multi-resolution segmentation. The original EM is applied on 
each of the original, parent and grandparent images. Each classification result is given a certain weight. Then, 
the final classification is calculated by combining the three weighted results as: 


Cay) = 0.4 * Covey) + 0.35 * City) + 0.25 * Cox) 


where, Co, Cı and C, are the results of the classification at pixel (x,y) for the original, parent and 
grandparent images respectively. Figure 2 shows the segmentation results of the WMEM, Figure 2(a) the 
original tumor image, Figure 2(b) the original segmented image using EM, Figure 2(c) the segmented parent 
image, Figure 2(d) the segmented grandparent image, Figure 2(e) the final segmented image using WMEM, 
and Figure 2(f) the extracted tumor. 


(d) (f) 


Figure 2. Stages of tumor extraction: (a) the original tumor image, (b) the original segmented image using 
EM, (c) the segmented parent image, (d) the segmented grandparent image, (e) the final segmented image 
using WMEM, and (f) the extracted tumor 


3.2. Tumor extraction 

Thresholding was used to separate the other segmented brain components (grey matter, white 
matter, and CSF) from the tumor. The thresholding produced images containing both tumor and skull 
because they almost have the same color intensity. So, the solidity feature was used to exclude the skull from 
the final image. Solidity is the proportion of the pixels in the convex hull that are also in the region. It is 
computed as Area/ConvexArea. It was found that tumors have higher solidity than the skull, so we selected 
the high solidity to get a binary image that only contains the tumor to be used in the reconstruction step [27]. 


3.3. Intermediate slice estimation 

Our aim is to visualize the brain tumor given a series of MRI brain slices. Our method is a 
modification of the classic shape-based interpolation technique [33] which is based on interpolating the 
values of the 3D volume that lay in the gap areas during slices acquisition. Shape-based interpolation is 
performed in two main steps. The first step is transforming the binary tumor image into a distance image. The 
second step is interpolating the distance image using a suitable interpolation function. 

Every series of MRI image is described with some important information such as the slice thickness, 
inter-slice distance and the in-plane resolution. These properties vary according to the MRI scanner strength 
and acquisition process. Before starting the interpolation, the number of gap slices should be calculated using 
the provided dataset properties. We computed the gap area in pixels by: 
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slice spacing 
in — plane resolution 


where, slice spacing is the gap between every two consecutive slices in millimeters and the in-plane 
resolution is the area of each pixel in millimeters [34]. 


3.4. Distance image transform 

The distance image represents the shortest distance from every pixel in the image to the nearest 
boundary pixel. Pixels inside the region of interest (tumor) are assigned with positive values, while pixels of 
the background are assigned with negative values. We calculated the distance transform by applying the 
following steps twice for each of the tumor and the background separately, then the two distance images are 
added to form the final distance image. 

Steps: 

Step 1: Create two distance images D and D'for the tumor and the background respectively, then set all the 
pixels of D and D’ to zero. 

Step 2: Create a counter I that increments with every iteration. 

Step 3: Starting from the boundary pixels, perform an erosion operation to the image to remove a single pixel 
from the tumor boundaries then add the value of J to the pixels of D at the locations of the eroded 
pixels. 

Step 4: Repeat step 3 until all the pixels of D are calculated. 

Step 5: Similarly, perform dilation to the tumor image by adding a single pixel starting from the boundary 
pixels and subtract the value of J in each iteration and insert its value in D’ at the locations of the 
dilated pixels. 

Step 6: Repeat step 5 until all the values of D’ are calculated. 

Step 7: By adding D to D’, the final distance image is formed. 

The previous steps are performed for every slice of the input volume. Figure 3 shows a sample of 
the distance image that results from applying the previous steps on one of our dataset images, the in-lined 
part represents the tumor area with positive values, while the outlined negative pixels represent the 
background and the zero pixels are for the tumor boundaries. 
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Figure 3. A sample of the resulting distance image 


3.5. Slice interpolation 

To reconstruct the final tumor volume, the gap area has to be interpolated using the given slices. 
Several interpolation methods have been used in medical imaging. Methods such as nearest neighbour and 
linear interpolation are fast and simple, yet they do not produce reliable results [35]. In our work, B-spline 
interpolation [36] was used on distance images as an input, and it proved to show accurate results. 

A separable spline model is chosen for interpolation. The p-dimensional spline function 
s(X),X = (x1,..,Xp) E RP is represented by the expansion 


S(X) = Xrezr c(k) p(X — k) (1) 
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where, c(k) are the B-spline coefficients. The basis functions in (1) are the integer shifts of the separable B- 
spline (X), which is a tensor product of the univariate B-splines of degree n: 


P(X) = P” Ca) B” Cp) (2) 


From the fold convolution of the box function, we can obtain the univariate B-splines. Their closed form 
expression is: 


Antin 


pra) = 3) 


where, x? = max{0,x}" is the one-sided power function and A"*? is the (n + 1) iteration of the central 
finite difference operator Af (x) = f (x + 2) =- f(x- 5. For a multidimensional image array f (i),i € ZP, 


the basic interpolation problem is to determine the coefficients c(k) in (1) such that the spline S (X) fits the 
pixel values exactly: S(X)|,-; = f(i),i E Z”. By applying this constraint and resampling (1) at the integers 
we get: 


FO = Xreze ck) eC — k). (4) 


Since (4) has the form of discrete convolution, we can determine the c(k) values by deconvolving the 
equation. In this study, experimental results showed that cubic spline interpolation proved to have higher 
accuracy than many interpolation methods such as linear, cubic and Hermite-splines. 


3.6. Interpolated slice binarization 

After applying interpolation on the input distance images of the tumor slices, the resulting 
interpolated distance images must be converted to binary images so that the gap slices can be inserted 
accurately. That was simply done be assigning all the zero and positive values to one (representing the tumor 
area) and zero to all the negative values (representing the background). By analyzing the results of 
shape-based interpolation, it was found that the interpolated image is always smaller than the ground truth 
image (i.e., the majority of error is produced as false negative pixels). This happens due to assigning all the 
interpolated negative values as background pixels, even if their value is very small. To overcome this 
drawback, we added a thresholding step to approximate the negative values so that they are not considered as 
background. The threshold was selected upon experimental trials. This modification improved the accuracy 
of the original method for all the used dataset images. 


3.7. Visualization 

After interpolating all the intermediate slices for the whole tumor images, the original and 
interpolated slices are stacked to form a 3D volume. Then, Isosurface is extracted and displayed in a 
3-dimensional view. Figure 4 shows the visualized brain tumor before slices interpolation, while Figure 5 
shows the visualized tumor using slices interpolation. It is noticed that the details of the tumor shape are 
displayed clearly after applying the proposed method. 


Figure 4. Displayed tumor volume without interpolation 
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Figure 5. The displayed interpolated tumor volume 


4. RESULTS AND DISCUSSION 

Data from the Internet brain segmentation repository (IBSR) was used in this study [37] which is 
provided by the Center for Morphometric Analysis at Massachusetts General Hospital. The data consists of 
multiple brain MRI scans for a patient with a tumor taken at roughly 6-month intervals over three and a half 
years. The algorithm was implemented using MATLAB® R2017b. Table 1 demonstrates the properties of 
the used MRI images series. 

As the provided data does not contain ground truth values, the system accuracy was calculated by 
removing an original slice from the MRI images, then interpolating it. The interpolated image was then 
compared with the original removed image and the percent error was computed as (5): 


|OA-IA| 
= — X 


6 100 (5) 


where, OA is the area of the object cross section in the removed slice and IA is the area of the object cross 
section in the interpolated slice. We compared our method with the original shape-based interpolation 
proposed by Raya and Udupa [33] and with Hermite interpolation used in [2]. 

The results are summarized in Tables 2 to 5. The results represent the accuracy percentage for each 
slice which is computed as 100 — 6, where 6 is the error percentage mentioned in (5). The first and last 
slices were excluded from the computations because for reliable interpolation results, the query points should 
lie in between the control points [33]. 


Table 1. Properties of the used MRI images series 


Series no. No. of slices with tumor Slice Spacing (mm) __In-plane Resolution 


126_10 17 3 1.015625 
126_13 16 2.5 1.015625 
126_21 17 25 0.937500 
126_26 17 0 0.937500 


Table 2. Accuracy percentage for each interpolated slice for series 126_10 of IBSR dataset 
Slice no. Original shape-based __ Hermite interp. Modified shape-based interp. 


2 83.8 88.2 92.1 
3 84.0 91.2 92.7 
4 85.8 86.6 91.6 
5 86.3 88.3 91.9 
6 90.0 92.4 95.0 
7 87.9 91.5 93.4 
8 90.8 92.2 94.6 
9 86.6 89.2 94.5 
10 88.1 91.1 93.9 
11 88.4 91.3 95.7 
12 88.2 89.6 94.8 
13 86.6 88.0 92.0 
14 85.5 90.4 93.6 
15 83.9 91.1 94.9 
16 83.5 92.7 95.2 
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Table 3. Accuracy percentage for each interpolated slice for series 126_13 of IBSR dataset 


Slice no. Original shape-based Hermite interp. Modified shape-based interp. 
2 19:2 81.4 86.7 
3 80.1 90.6 91.6 
4 81.3 88.0 92.3 
5 86.5 92.1 94.1 
6 87.0 92.1 94.3 
7 88.7 92.1 94.6 
8 88.8 94.9 96.0 
9 87.7 91.0 95.8 
10 88.7 95.0 96.1 
11 88.6 90.7 94.7 
12 89.7 94.8 96.8 
13 87.6 89.8 92.9 
14 81.4 83.2 90.0 
15 79.1 84.6 89.9 
mean +variance 85.31+16.53 90.02+18.50 93.27+8.382 


Table 4. Accuracy percentage for each interpolated slice for series 126_21 of IBSR dataset 


Slice no. Original shape-based_ Hermite interp. Modified shape-based interp. 
2 83.3 87.6 92.4 
3 83.0 85.5 91.6 
4 84.7 87.6 92.1 
5 86.3 91.9 92.6 
6 92.3 94.5 96.0 
y 92.8 96.2 97.1 
8 90.9 93.9 96.2 
9 92.4 93.0 95.9 
10 89.6 90.5 93.8 
11 90.0 89.5 92.8 
12 91.7 95.6 97.4 
13 88.6 88.9 95.9 
14 92.1 95.9 97.0 
15 89.2 91.5 93.9 
16 81.7 85.8 89.7 
mean +variance 88.57+14.556 91.19+13.19 94.29+5.629 


Table 5 Accuracy percentage for each interpolated slice for series 126_26 of IBSR dataset 


Slice no. Original shape-based Hermite interp. Modified shape-based interp. 
2 83.4 88.9 92.2 
3 89.4 90.7 93.9 
4 90.4 91.4 94.5 
5 91.3 92.0 95.6 
6 89.8 92.3 95.8 
7 89.6 93.7 95.6 
8 87.9 91.4 94.9 
9 92.9 94.2 97.1 
10 93.0 95.2 97.7 
11 92.2 93.8 96.7 
12 92.8 93.6 96.1 
13 90.5 93.4 95.8 
14 91.4 93.2 95.9 
15 90.6 93.1 96.3 
16 85.7 89.9 92.7 
mean +variance 90.06+7.248 92.45+2.987 95.38+2.33 


The final accuracy of the interpolated tumor volume is expected to be higher than the accuracy of 
the previous results. As when removing an original slice, the distance between the given slices (i.e., control 
points) is doubled, which will eventually decrease the interpolation accuracy. It was also noticed in Figure 6 
that the error pixels are distributed evenly around the tumor boundary, Figure 6(a) the original tumor image, 
Figure 6(b) the interpolated tumor image, and Figure 6(c) the difference between A and B, green pixels 
represent the false positives and red pixels represent the false negatives. We used the structural similarity 
index (SSIM) to compare the interpolated and the original images. SSIM tells how much one image is 
structurally similar to the other [38], [39]. Values between 98% and 99% were gained. This shows that the 
interpolated tumor preserved its original shape and curvature areas, which will not affect the medical analysis 
of the system output. Since no automatic method has yet been proved to be reliable for usage as ground truth, 
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human expertise is still needed to confirm the efficiency of the automatic methods [40]. We contacted some 
experts in the Neurosurgical department of Ain Shams University, and they approved the reliability and 
efficiency of the reconstructed tumor images. 


D 


(a) (b) (c) 


Figure 6. Difference between original and interpolated tumor images, (a) the original tumor image, 
(b) the interpolated tumor image, and (c) the difference between a and b, green pixels represent the false 
positives and red pixels represent the false negatives 


5. CONCLUSION 

Brain tumor visualization can efficiently help in treatment planning and resection assessment of 
brain cancer. Having accurate information about the shape and size of the tumor will increase the chances of 
total resection, which will eventually increase the mean survival period of the patient. Many attempts have 
been done for brain tumor visualization using different methods such as alpha-shapes, marching cubes, point- 
based methods, and interpolation-based methods. 

In this paper, we presented a modified shape-based interpolation method for 3D reconstruction of 
brain tumors using a group of parallel MRI image slices. The tumor was first segmented from each brain 
slice using the WMEM algorithm, then the gap area was computed using our method. The final 3D shape of 
tumor was then displayed using isosurface function. The accuracy of the system was calculated by removing 
an original slice from the series and interpolating it, then comparing both images. An average accuracy of 
94.1% was gained when comparing the images pixel by pixel, while the average SSIM was 98.5% which 
reveals that the proposed method preserves the structural shape of the interpolated slices. 
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